Fluctuation formula and non-Gaussian distribution in the thermostated Lorentz gas 
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In this paper we numerically examine the connection of the Gallavotti-Cohen fluctuation formula 
and the functional form of the corresponding probability density function in the field driven Lorentz 
gas thermostated by the Gaussian isokinetic thermostat. We analyze the moments of the entropy 
production rate fluctuations and show that all the central moments of the averaged fluctuations 
exhibit power law dependence on the length of the averaging time interval, indicating that this 
density deviates from a Gaussian. Furthermore the obtained exponents are found to obey a special 
pairing rule showing that the corresponding probability density function can not be scaled in the 
averaging time interval. 
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I. INTRODUCTION 

In recent years a large number of papers have focused 
on the Gallavotti-Cohen fluctuation formula (FF) with 
both theoretical and numerical tools in various nonequi- 
librium systems B S S S S B S S H E3 The FF is 
a symmetry property of the probability density function 
(PDF) of a dynamically measured quantity £ connect- 
ing the probabilities of measuring £ values with equal 
magnitudes but opposite signs. This property was first 
observed numerically in a system of thermostated fluid 
particles undergoing shear flow Q] where £ was the phase 
space contraction rate. The significance of this result was 
largely affected by the fact that the phase space contrac- 
tion rate could be connected to the entropy production 
rate, thus the FF gave an insight to the noncquilibrium 
thcrmodynamical behavior of the system. Since then the 
nontrivial connection of the phase space contraction rate 
and the entropy production rate has been investigated 
and was treated in [Tlj . 

More precisely, let £ T (t) denote the quantity £ averaged 
over a time interval of length r centered around time t: 

r/2 

£ T (t) = - J £(t + t')dt'. Considering it as a stochastic 

T -r/2 

variable, its statistical properties in a steady state can 
be characterized by a probability density function S r (x) . 
The FF states that this PDF has the following property: 
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In other words, the probability of observing a negative 
£ value is exponentially smaller than that of the corre- 
sponding positive value in the r — > oo limit. 

The FF has also been found to be valid in certain sys- 
tems with strong external forcing ^1 > thus it is expected 
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to say something important about systems far from equi- 
librium (as opposed to linear response theory, which is 
valid for vanishing external fields). Consequently, signifi- 
cant theoretical efforts have been made to find a common 
propertybehind this behavior in all the systems obeying 
the FF BUSH 

Since the FF connects the values of the S T (x) func- 
tion on the negative half axis to the values of S T (x) on 
the positive one, but it does not state anything about 
the functional form of S T (x) , it is interesting to examine 
the latter in systems found to obey the FF. This ques- 
tion is emphasized by the observation that Eq. Q can 
be satisfied by a Gaussian PDF with a subcondition con- 
necting its average and variance. This raises the necessity 
to study how close the underlying PDF of S T (x) is to a 
Gaussian in systems which seem to obey the FF. 

In a previous paper we investigated numerically 
the validity of the FF in the two and three dimensional 
periodic LG subject to a constant electric and magnetic 
field and thermostated by the GIK thermostat. We found 
that the symmetry property described by the FF is valid 
for low and high r values in the reversible configurations 
and seems to be approximately valid in the r — > oo limit 
even in the absence of reversibility. This also motivated 
us to study the r-dependence of the details of S T (x) in 
the thermostated LG. 

In this paper we numerically analyze the connection 
of obeying the FF and the details of the functional form 
of the corresponding PDF in the periodic Lorentz Gas 
(LG) thermostated by a Gaussian Isokinetic (GIK) ther- 
mostat. In Sec. [n] we describe the simulated system, in 
Sec. IIIII we present our numerical results and in Sec. II VI 
we summarize our conclusions. 



II. THE SYSTEM 

One of the simplest and most investigated models suit- 
able for studying transport phenomena is the field driven 
Lorentz gas (LG) thermostated by a Gaussian isokinetic 
(GIK) thermostat. This model consists of a charged par- 
ticle subjected to an electric field moving in the lattice of 
elastic scatterers. Due to the applied electric field, one 
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must use a thermostating mechanism to achieve a steady 
state in the system. Such a tool is the Gaussian isoki- 
netic thermostat which preserves the kinetic energy of 
the particle [l4j . 

For the sake of simplicity, we investigate in our study 
the two dimensional (2D) periodic Lorentz gas with 
circular scatterers arranged on a square lattice. We 
present the equations of motion in dimcnsionless vari- 
ables: choosing the unit of mass m and electric charge q 
to be equal to the mass and electric charge of the parti- 
cle, yields m = q = 1 in our model. The unit of distance 
is taken to be equal to the radius of scatterers (R = 1), 
and the unit of time is chosen to normalize the magni- 
tude of particle velocity to unity. Let q = (<?i, 92) denote 
the position and p = (pi,P2) the momentum of the par- 
ticle (|p| = 1). The phase space variable of the system 
is r = (q, p) ; it is transformed abruptly at every elastic 
collision and evolved smoothly by the differential equa- 
tion 

q = P 

p = E-ap (2) 

between the collisions. Here a is called the thermostat 
variable, while E is the constant vector playing the role of 
the external electric field. The GIK thermostat requires 
the fixing of the kinetic energy, which leads to the choice 
a = Ep in Eq. J2J. 

Dissipation in such models can be measured by the 
phase space contraction rate a. It can be computed by 
taking the divergence of the right-hand side of Eq. J3J : 

cr = -divf = a, (3) 
and it can be shown that in this model 

a(t) = £(t), (4) 
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FIG. 1: The probability density II T (a;) in the 2D periodic 
LG, with lattice constant d = 2.1 and external field E = 
(0.3,0.48), the number of collisions is 1.58 x 10 7 and the av- 
eraged time between two collisions is ~ 0.62. For higher r 
values the curves seems to be "Gaussian like" : fitting a Gaus- 
sian function onto the measured values yields very good visual 
agreement. 

of possible £ r values depends on the strength of the elec- 
tric field. Instead of £ T , we may introduce the quantity 

2 

n T (t) = i f n E p(t + t')dt, (5) 

where denotes the unit vector parallel to E. Since the 
magnitude of p is unity, tt t always satisfies ir T e [—1,1]. 
By making a histogram of the periodically measured val- 
ues of tt t , one gets an approximation of its probability 
density II T (x) shown in Fig.^ The two stochastic vari- 
ables defined above satisfy £ r = Ett t , so the relationship 
between two probability densities S T and II r is simply 

s T (z) = -|n T (t). 



where £ is the entropy production rate defined accord- 
ing to the £ = expression of irreversible thermody- 
namics, and the temperature T is given by analogy with 
kinetic theory [ll|. We note that while Eq. Q is valid 
for the GIK thermostat, it need not be true in general, 
as shown e.g. in models thermostated by some variants 
of the Nose-Hoover thermostat |l£| or by deterministic 
scatterings [Trij . 



III. NUMERICAL RESULTS 

The goal of our numerical simulations is to measure 
H x (x) with a precision that is sufficient to check the va- 
lidity of the fluctuation formula, as well as the question 
of Gaussianity of the distribution. Due to Eq. S T (x) 
can be measured by periodically computing £ T along a 
long particle trajectory and making a histogram of these 
data. The disadvantage of this method is that the range 



A. Scaling 

We can observe in Fig.^that for high enough r values 
the PDF n T (a;) seems to be close to a Gaussian. In order 
to examine this behavior we applied the following scaling 
transformation to the II T (a;) functions: 



S T (x) 



1 



(max) 



n 



(G) 



In this equation x max denotes the x- value where II r (x) 



takes its maximal value denoted by Yr r 



(max) 



Observing 



in Fig. |3| that the standard deviation of n T (x) is propor- 
tional to \j\fr (thus Y[[™- ax% > y/r), one can see that 
Fig. El shows 7iv fluctuations of 0(1/ y/r). It is worth 
noting that fluctuations of this magnitude also constitute 
the range where the central limit theorem (CLT) could in 
principle guarantee the Gaussian form of S T (x) for high 
enough r values. Indeed, by examining Fig. [21 we can 
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FIG. 2: The S T (x) scaled function of IF (a;) for different r 
values with field E = (0.3, 0.48) (the same configuration as 
for Fig. . The dashed curve shows the Gaussian function 
S G {x) = e-^ 2 . 

conclude that for high r values all the measured points 
of S T (x) are close to the Gaussian Sg{x) — e~^** x2 . 

B. The moments 

In order to characterize the deviation of the PDF from 
a Gaussian, we compute its moments and cumulants from 
the simulated distribution. Let m™, c'*) and denote 
the i'th order moment, central moment and cumulant of 
II T (x) respectively (i = 1,2,3,...) [l7l| . First we deal 
with the central moments of II r (x) for a given configura- 
tion of the periodic LG, shown in Fig. [3J Based on this 
figure, we can make several important observations. 

1. The central moments of order i exhibit power 
law dependence on r with exponents /3, depending 
on i: Cr ~ c^r - ^*; Table (JU lists the numerical 
values of ai and Pi found in the simulation. 

2. The exponents /3, are found to form pairs, at least 
approximately: p^ ~ (3^, P& ~ Pq, and so on. 

3. The even order central moments behave the same 
way as for a Gaussian distribution, i.e., they all can 
be computed from the variance alone. The Gaus- 
sian moments Cq calculated in this way from c^ 2 -* 
are shown as continuous lines in the figure, and 
indeed, they seem to run very close to the actual 
measured values of the even order central moments. 

Analyzing several configurations with various strengths 
and directions of the external field, we found that the Pi 
exponents are approximately the same in all the config- 
urations. In contrast, the ai factors do depend on the 
details of the chosen configuration. 

Now we compare our findings to the expectations based 
on a Gaussianity assumption. Because of the power law 
dependence in the odd order central moments, we can say 



FIG. 3: The central moments of IF (a;) in the 2D periodic 
LG. In this configuration E = (0.3, 0.48), the collision number 
is 1.5 x 10 7 . The continuous lines show the behavior expected 
from Gaussian moments: c^q = 3 c' 2 ' , eff = 15 c' 2 - 1 , = 
105c (2)4 , 4 10) = 945c (2)5 , etc., calculated from c (2) . 
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4.937 ±0.038 f= 


a 5 


3.742 ± 0.656 


10 


5.068 ±0.187 f= 


a 5 


2.631 ±0.187 



TABLE I: The exponent j3i and the prefactor ai of the i-th 
central moment shown on Fig. These values are computed 
using the NLLS Marquardt-Levenberg algorithm included in 
gnuplot to fit the measured data points with a power function 
on the interval [100, 1000]. The errors presented here come 
from the fitting method. The odd order moments seem to 
run close to their even successors: fan+i ~ 02(n+i) ~ n + 1- 



that n r (x) is definitely not Gaussian for any finite value 
of r (they all should be zero for a Gaussian distribution). 
On the other hand, the distributions n T (x) do seem to 
converge to a Gaussian in the limit r — > oo. This can be 
checked in Fig.0|showing the behavior of the cumulants: 
all of the higher (i > 2) order cumulants go to zero faster 
than those of order i = 1, 2 in the t-ioo limit. 

Another important consequence of the pairing of the 
exponents is that H T (x) can not be scaled in the sense of 
Eq. 0. Indeed, the equation S T (x) = S(x) results that 

LI r (x ± m) = Ii i f nax ^ S (rii max ^ x^j , since x max = m for 
high r values. Substituting this equation into the defini- 
tion of the central moments and performing the variable 
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FIG. 4: The absolute values of the cumulants of II T (a;) in 
the periodic LG. The configuration is the same as for Fig. |3 
H T (x) converges to a Gaussian in the r — > oo limit. Since 
g( 3 ) r*j r~ 2 , the convergence of IT- (a;) to a Gaussian is power 
like; the deviation is on the order of r -2 . 



transformation y = U^^x yields 



/oo 
dx!l T (x + m)x l = 
-OO 



The last expression shows that if S T (x) = S(x) were true 
in this system, then the /3, exponents should obey the 
(3{ = 4 identity, since ni™ ^ ~ ^/r for high r values. 
Table Q] shows that the odd order /3j exponents deviate 
significantly from |. This behavior indicates that the 
scaling property suggested by Fig. can only be treated 
as a good approximation and not as an exact result in 
the GIK thermostated LG. This also raises the question 
whether similar claims of scaling, based on the visual 
properties of the PDF in other systems as e.g. in Ref. [Isj . 
could also be supported by the behavior of the moments. 



The Fluctuation Formula 



We can rewrite the FF in the terms of II T (x) as 



lim II la=5£L 

T^oo E T n r ( — X) 



(8) 



In order to visualize the FF, we introduce the quantity 

(9) 



11, II T (a;) 
D T (x) = — - In- v ; 



Et U t (-x) 

If the FF is satisfied in the system, then D T (x) is linear 
with slope 1 in the r — + oo limit. Analyzing the behavior 
of this quantity in Fig. one can make two important 
observations: 
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FIG. 5: D T (x) in the same configuration as for Fig. The 
inset shows D T (x) for higher r values. The inset and the 
figure have the same axes. 



1. The FF is satisfied quite well for both low t val- 
ues and for large ones (in this context, low t value 
means 6-8 times the average time between two sub- 
sequent collisions). 

2. As r gets higher, the interval where D T (x) can be 
plotted is getting narrower and narrower. This is 
a consequence of the fact that for higher r values 
large fluctuations of £ T become less and less proba- 
ble. In other words fluctuations outside the plotted 
range are practically unobservable. 

We note that we tested several other configurations with 
various choices of the magnitude and direction of the 
electric field E in the ergodic range and we didn't 
find any relevant deviation from the behavior presented 
above. 



D. The Gaussian approximation 

Since we have seen that the even order central mo- 
ments behave in a Gaussian way, it is worth checking 
if the subcondition that leads to the FF for Gaussian 
distributions holds for LI T (a;). This subcondition can be 



derived by substituting II T (x) = — exp ■ — 
into Eq. JHJ. This leads to the following formula: 



V 2 -? ) 



11 2m T 
hm 5- 



= 1 . 



(10) 



where m T = and o\ — S^ 1 according to the previous 
notation. 

In Fig.Elwe can verify that the jj-f^pr = 1 identity is 
violated for all finite r values but satisfied in the r — > 00 
limit. On the other hand, we have seen in Fig. that the 
FF is valid with reasonable precision even for low values 
of r. These observations indicate that: 

1. For finite r values the the odd order moments also 
contribute to the fulfillment of the FF, thus the 
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FIG. 6: The quantity at different field strengths E 

in the 2D periodic Lorentz Gas. In these configurations E is 
parallel with (0.3,0.48) and the collision number is « 1.5xl0 7 . 
The curves deviate from 1 for all the r values presented on 
the figure, however seem to converge to 1 in the r — > 00 limit. 



E. The observability of large fluctuations 

Since we have shown that for high r values neither 
n T (a;) nor D T (x) can be observed on the whole inter- 
val [—1,1], we introduce A t ^e to quantify where D T (x) 
can be observed, i.e. where the FF can be tested. Let 
[— A T: e, A T ,s] denote the interval where D T {x) can be 
reconstructed from a given time series in a specific con- 
figuration. Data from our numerical simulations are sum- 
marized in Fig. It shows that A t< e has a power-law 
dependence on r: 

A r , B ~ r- 5 ^ , (11) 

where the exponent <5 tends to increase with increasing 
field strengths. 



Gaussian approximation of the PDF Il r (a;) violates 
the FF. This behavior is also supported by the anal- 
ysis performed in Ref. 6J. 

2. For high r values the observable part of the PDF 
Il T (a;) seems to converge to a Gaussian (see. Fig- El) 
and the Gaussian approximation of the PDF H T (x) 
satisfies the FF. 

This is in agreement with the observed behavior of the 
cumulants indicating a convergence to Gaussian distribu- 
tions for r — * 00; indeed, for Gaussian distributions the 
FF is equivalent to Eq. (|rU|) . 




FIG. 7: The dependence of A T: e on r for different field 
strengths. The direction of E is parallel with (5, 8) but its 
magnitude varies. The curves appear to be linear in the dom- 
inant region on the log-log plot, suggesting a power law de- 
pendence on r. 



IV. CONCLUSION 

In this paper we have investigated numerically the 
PDF of the EPR fluctuations for both low and high aver- 
aging time intervals in the 2D periodic LG thermostated 
by the GIK thermostat. We have shown that for low av- 
eraging time intervals the EPR fluctuations are definitely 
not Gaussian, but they satisfy the Gallavotti-Cohen fluc- 
tuation formula. For high averaging time intervals the 
EPR fluctuations are converging to a Gaussian distribu- 
tion and still obeying the FF, thus the Gaussian approx- 
imation of the EPR fluctuations seems to be better and 
better as t gets higher in this model. We have analyzed 
the moments of the corresponding PDF to characterize 
its functional form and found a special power-law depen- 
dence (pairing of the (3i exponents) of the central mo- 
ments on the averaging time interval. Our results raise 
the question whether the /^-pairing could be observed in 
other systems satisfying the FF, or this is just a special 
property of the GIK thermostated Lorentz gas. 



Acknowledgments 

The authors are grateful to Tamas Tel for fruitful dis- 
cussions and a careful reading of the manuscript. Fur- 
thermore the first author thanks Jurgen Vollmer for valu- 
able e-mail conversations. This work was supported 
by the Bolyai Janos Research Grant of the Hungarian 
Academy of Sciences and by the Hungarian Scientific Re- 
search Foundation (Grant No. OTKA T032981). 



[1] D.J. Evans, E.G.D. Cohen and G.P. Morriss, Phys. Rev. [2] G. Gallavotti and E.G.D. Cohen, Phys. Rev. Lett. 74, 
Lett. 71, 2401 (1993) 2694 (1995) 



6 



[3] G. Gallavotti, Phys. Rev. Lett. 77, 4334 (1996) 

[4] J. Kurchan, J. Phys. A, 31, 3719 (1998) 

[5] D.J. Searles and D.J. Evans, Phys. Rev. E 60, 159 (1999) 

[6] F. Bonetto, G. Gallavotti and P.L. Garrido, Physica 

105D, 226 (1997) 
[7] D.J. Evans and D.J. Searles, Phys. Rev. E 52, 5839 

(1995) 

[8] C. Maes, J. Stat. Phys, 95, 367 (1999) 
[9] L. Rondoni and E. Segre, preprint chao-dyn/9810028 
[10] L. Rondoni, T. Tel and J. Vollmer, Phys. Rev. E 61, 
4679 (2000) 

[11] See e.g. E.G.D. Cohen and L. Rondoni, Chaos 8, 357 
(1998) 

[12] For an overview, see G. Gallavotti, Chaos 8, 384 (1998) 

and references therein. 
[13] M. Dolowschiak and Z. Kovacs, Phys. Rev. E 66, 066217 

(2002) 

[14] For a review of deterministic thermostats, see G.P. Mor- 
riss and CP. Dettmann, Chaos 8, 321 (1998) and refer- 
ences therein. 

[15] K. Rateitschak, R. Klages and W.G. Hoover, J. Stat. 
Phys. 101, 61-77 (2000) 



[16] K. Rateitschak and R. Klages, Phys. Rev. E 65, 036209 
(2002) 

[17] Obviously the quantity in question is H r (x), however it is 
easier to work with Il r (a;). The connection between the 
analogous quantities of the two PDFs S T (x) and II T (a;) 
can easily be derived. Let rh^ % \ and q^ denote the 
i'th order moment, central moment and cumulant per- 
taining to E T (x), then: m (i) = E l m (i) , c (i) = ETc^ and 
= E l q^. 

[18] L. Rondoni and G.P. Morriss, Preprint ESI 1148 (Vienna, 
2002) 

[19] The dynamics of the periodic LG thermostated by a 
GIK thermostat is nonergodic for high external fields, as 
shown first in J.Lloyd, L. Rondoni and G.P.Morriss, Phys. 
Rev. E 50, 3416 (1983). Since we have used a different 
lattice (cubic) in our numerical study, we checked the er- 
godicity of our model for different field values and found 
behaviors similar to the ones reported in that paper: for 
low E values the dynamics seems to be ergodic, but for 
high E values (E > 6) the breakdown of ergodicity can 
be observed. 



